Insights into planktonic food-web dynamics through the lens of size and season

Knowledge of the trophic structure and variability of planktonic communities is a key factor in understanding food-web dynamics and energy transfer from zooplankton to higher trophic levels. In this study, we investigated how stable isotopes of mesozooplankton species varied seasonally (winter, spring, autumn) in relation to environmental factors and plankton size classes in a temperate coastal ecosystem. Our results showed that spring is characterized by the strongest vertical and size-structured plankton food-web, mainly fueled by the phytoplankton bloom. As a result, spring displayed the largest isotopic niche space and trophic divergence among species. On the contrary, both pelagic and benthic-derived carbon influenced low productive seasons (winter and autumn), resulting in more generalist strategies (trophic redundancy). Stable isotope mixing models were used to explore how different seasonal structures influenced the overall food web up to predatory plankton (i.e., mysids, chaetognaths, and fish larvae). Different feeding strategies were found in spring, with predators having either a clear preference for larger prey items (> 1 mm, for herring and dab larvae) or a more generalist diet (sprat and dragonets larvae). During low productive seasons, predators seemed to be more opportunistic, feeding on a wide range of size classes but focusing on smaller prey. Overall, the food-web architecture of plankton displayed different seasonal patterns linked to components at the base of the food web that shaped the main energy fluxes, either from phytoplankton or recycled material. Additionally, these patterns extended to carnivorous plankton, such as fish larvae, emphasizing the importance of bottom-up processes.


Study area
The Eastern English Channel and Southern Bight of the North are epicontinental seas bordered by France, Belgium, and England (Fig. 1) subject to high anthropogenic pressures (marine traffic, fishing activities, gravel extractions 33,34 ).This area is very shallow, and highly productive, with strong vertical and horizontal mixing, and a seasonal temperature gradient 35 The area has been extensively studied for more than a century to advise management decisions relating to environmental health, sustainable use of resources, and conservation 36 .With the objective of an ecosystem-based management, a wide range of available data on both abiotic (e.g., tidal hydrodynamics, sediments, Chl a) and biotic factors (e.g., species abundance, distribution and composition) have been used to implement ecosystem-based models in the EEC 37,38 .This holistic view of the ecosystem requires knowledge of species interactions in terms of predator-prey relationships together with food-web structure and functioning.For the most part, these interactions are well documented and informed for commercial fish species, but a lack of understanding remains for lower trophic levels 12,13 .Bottom-up effects, however, have recently been identified in the EEC as one of the most important drivers of variation for outputs in ecosystem-based models, affecting ecosystem dynamics and changes in biomass across all functional groups and trophic levels 39 .Similarly, studies on the community composition in the North Sea have identified that zooplankton composition and abundance are main bottom-up drivers of ecosystem dynamics 40,41 .In the EEC, spatial data comes primarily from fisheries-oriented oceanographic surveys that usually neglect the zooplankton compartment, or from dedicated ichthyoplankton surveys (e.g., International herring larvae survey, IHLS), which are restricted in time to just one season 42 .Seasonal dynamics of zooplankton remain largely uninvestigated in the EEC, with data only available from the French coastal station of Gravelines [43][44][45] .

Sampling procedures
Samples were collected opportunistically from 2017 to 2019 in winter (January to February from the International Bottom Trawl Surveys, https://doi.org/10.18142/17),spring (March to May 2017, from the REIVE I and II surveys, https:// doi.org/ 10. 17600/ 17010 400, and PHYCO surveys https:// doi.org/ 10. 17600/ 17010 500) and autumn (September to October 2017 to 2019, from the Channel Ground Fish Surveys https://doi.org/10.18142/11)(Fig. 1).Environmental data for each survey and additional details on sampling protocols can be accessed through the survey's DOI.Briefly, environmental parameters displayed a strong seasonal pattern with temperatures ranging from a minimum of 6.5 °C in winter to a maximum of 18 °C in autumn.Salinity was rather stable throughout the seasons (~ 34) with minimum values found at stations close to estuary mouths (the Seine, Somme, Authie, and Canche estuaries).Correspondingly, the influence of plumes of turbidity in front of estuaries in spring resulted in elevated average values of suspended particular matter (SPM) during this season with high variability.As expected, values of Chl a were higher in spring than in winter and autumn.Conversely, dissolved nutrients (NO 2 , NO 3 , PO 4 , SiOH 4 ) were lower in spring compared with the other seasons (Supplementary Table S1).Niskin bottles were used to collect surface (1 m depth) water samples.Samples were immediately pre-filtered on a nylon-mesh filter of 200 µm to remove large zooplankton.The remaining fraction (< 200 µm) was then filtered through pre-combusted (450 °C for 4 h) GF/F filters until saturation (usually 1-2 L).Mesozooplankton was collected from vertical hauls (from 4 to 58 m depth) using a WPII net (200 µm mesh size).Fish larvae were caught using either a midwater ring net (winter and autumn) or a bongo net (spring) both with a mesh size of 500 µm.Fish larvae were sorted on board.All samples were immediately frozen at − 80 °C.All methods were carried out in accordance with relevant guidelines and regulation.
Stable isotope compositions (δ 13 C, δ 15 N values) from the main plankton species were used to decipher the primary energy pathways.In the laboratory, zooplankton and fish larvae were rinsed in distilled water, sorted, measured and identified to the lowest possible taxonomic level 46,47 .One to ~ 100 individuals of each taxon of similar size-classes were pooled together to ensure sufficient mass (~ 300 µg) for stable isotope analysis.Five to ten random individual measurements were taken from each pool to obtain average ranges of total lengths (mm).Samples were then freeze-dried, homogenized and ground to a fine powder.GF/F filters were observed under a stereomicroscope to remove nauplii and small zooplankton if present.Every filter was split in half and carbonates were removed from one half by fuming with HCl for subsequent δ 13 C analysis, the other half was kept for δ 15 N analysis.Filters were then freeze-dried to remove any excess water.Isotope ratios were measured with a Thermo Delta V isotope mass ratio spectrometer, coupled to a Carlo Erba NC 2500 elemental analyzer.The accuracy of the isotopic ratio measurements was checked by repeated analyses of an in-house standard (one analysis of the standard after every 10 samples) with an overall standard deviation of 0.2‰ for both elements.Stable isotopes ratios were expressed following the classical δ notation with: where X being δ 13 C or δ 15 N, and R the isotopic ratios (13C/12C or 15N/14N, respectively) measured in samples and in international standards (Vienna Pee Dee Belemnite for C and atmospheric nitrogen for N).In total, the sampling effort resulted in 552 measurements: 289 from zooplankton (18 taxa), 167 from fish larvae (13 taxa) and 96 from GF/F filters of water samples (Table 1).

Data analysis
Sources of variations in stable isotopes and size-structure models One of the challenges when comparing food-webs over different spatial and temporal scales is the variation in baseline isotopes that needs to be accounted for when comparing consumers' values.For this, we investigated the relationship between δ 13 C and δ 15 N isotopes and spatial and temporal factors, using a Generalized Additive Model (GAM) with a spatial smooth term.The GAM is a flexible statistical approach that allows for the modeling of non-linear relationships, making it well-suited for capturing complex patterns in environmental data.Prior to model fitting, we transformed the δ 13 C values to be positive (δ 13 C positive ) by adding a constant value and then applied a log transformation (δ 13 C log ) to meet the normality assumptions required for the analysis.The GAM models included a tensor product smooth of the "longitude" and "latitude" variables for n = 96 seston samples to account for potential spatial autocorrelation, and a categorical variable for « Season» with three levels (i.e.winter, spring or autumn) to include any seasonal effects on baseline stable isotopes values.Uneven seasonal sampling, with spring limited to 2017 and autumn/winter spanning multiple years (2018, 2019), hindered the complete capture of yearly influences on isotopic baselines.Before modeling, the "year" impact on winter and autumn baselines was explored.However, statistically non-significant effects prompted the pooling of data across years for a holistic model representation.The models were fitted using the gam function in the mgcv package in R with the Gaussian family and REML method for estimating the model parameters.
where s(longitude,latitude, bs =" tp", k = 10) specifies a smooth term for the spatial coordinates, longitude and latitude, using a thin plate spline basis function with k degrees of freedom.The choice of k in our GAMs was pivotal to balance model complexity and effectively capture spatial variations in our isotopic data.In the carbon model, convergence challenges led to increasing knots to k = 30, addressing the intricacies in carbon isotope spatial variability without compromising model stability.Meanwhile, the nitrogen model used k = 10, striking a balance between capturing spatial patterns and maintaining model simplicity.The smooth term allows for flexible modeling of the spatial variation in the response variable.Convergence of the GAM models was assessed using several diagnostic plots (see Supplementary S1 for details).The predicted values, associated with particular coordinates and seasons, served as the isotopic baselines to adjust the observed isotopic data in zooplankton.This adjustment involved subtracting the predicted baselines from the observed zooplankton values.The resulting dataset, termed δ 13 C adjusted or δ 15 N adjusted , reflects the refined values based on the derived isotopic baselines.Due to the log transformation in the carbon model, a back-transformation was required before this correction was applied for δ 13 C.Further details on spatial patterns on baseline estimates can be found in the Supplementary S1.
The three seasons (winter, spring and autumn) can then be compared and represented on a two-dimensional plane, the two axes of which being the adjusted isotopic signatures of nitrogen (δ 15 N adjusted ) and carbon (δ 13 C adjusted ).For each plot, a convex hull, illustrates the overall theoretical niche space occupied by the plankton and is the equivalent of the richness isotopic functional diversity metric proposed by Cucherousset and Villéger  (2015).Details on formulas and calculations can be found in references 48,49 .The inside color polygons correspond to the niche space occupied by the plankton at each season.A larger seasonal polygon reflects a larger diversity of food resources and feeding strategies.At the species level, species close to the center of the seasonal polygon reflect more generalist species, and species at the edges reveal higher trophic divergence (i.e.specialized feeding preferences).
Linear mixed effect models (LMEM) were chosen to explore sources of variations (i.e.size, species, season) of stable isotope values for consumers.LMEM were particularly appropriate for the structure of our data, which encompass different grouping factors, unbalanced configuration with different sample sizes, and nested data (not truly independent data).All models were fit using the 'lmer' function in the "lme4" 50,51 and "vegan 52 packages in R version 4.1.0 51, with restricted maximum likelihood (REML) estimation used to estimate the model parameters.To ensure model convergence and validate assumptions of normality and homoscedasticity of residuals, diagnostic plots of the fitted model were examined by using the package 'performance' 53 .Likelihood ratio tests via anova() were used to test the significance of both, fixed factors (i.e.season, size) and random effect structures by comparing full models against reduced models (see Supplementary S2 for details).Full models included the interaction between "Season" and "Size" (log transformed) as a fixed factors, with "(log_size | Species)" as a random effect thus accounting for potential interspecific variations within the size effect.
Best models for carbon and nitrogen were: www.nature.com/scientificreports/

Energy fluxes
Food-web topology A food-web topology was defined for the EEC from phytoplankton to fish larvae.Following the general methodology proposed by Planque et al. ( 2014) 54 , the topology consisted of three elements: nodes (i.e.species), links (i.e.trophic interactions), and directions (i.e. who is the predator and who is the prey).Trophic interactions and directions were constructed based on data from either peer-reviewed publications, gray literature or institutional reports (documented interactions, TPlank0).These were completed by inference, on the basis of knowledge on similar species from comparable regions and maximum prey/predator ratios (length of the largest prey divided by the length of the predator 55 (potential interactions, TPlank1) (see Supplementary S3 for details).The presence of some groups in the plankton e.g.meroplankton, is limited to specific seasons or life stages on the basis of their specific traits, like fast growth or short life span.As a result, all species listed in the general topology do not necessarily meet and interact.Seasonal variations in plankton assemblages were investigated by constructing food-web topologies characteristic of winter, spring and autumn communities.

Stable isotope mixing models
Isotope mixing models are based on the principle that a consumer's isotope values result from the mixing of the isotope values of its food sources proportionally to their relative contributions to its diet (after adjustment for isotopic fractionation during digestion, metabolism, and assimilation, i.e. trophic enrichment factor TEF 56 ).In this study, we used a combination of two mixing models (IsoWeb 57 and MixSIAR 58 ) to identify drivers of variation in energy pathways for main consumers as proposed by Giraldo et al. ( 2017) 59 (see Supplementary S4 and S5 for details).MixSIAR is a consumer-scale mixing model and was used to explore variation (season, plankton size) in trophic pathways for predatory plankton.Prey sources were only considered if at least three samples were available.Permutational multivariate analysis of variance (PERMANOVA, 999 permutations) 60 was used to test differences in centroids and dispersion (based on δ 13 C and δ 15 N values) among prey sources using the adonis function in the vegan package in R. Prior to permanova, the homogeneity of dispersion among the different species was tested using the betadisper function.Multilevel pairwise comparison (with Bonferroni-corrected p values) was then used to identify when two prey sources were indistinguishable based on their isotopic signatures (pairwiseAdonis package) 61 .Sources were aggregated a priori only if (i) there were no differences on their centroid position, (ii) the combined sources had some functional or ecological significance, and (iii), they were of similar sizes.Following suggestions by 62 , visual inspection of the final isotopic space (i.e.predator values should fall within the isotopic space created by prey values after TEF corrections) and correlation coefficients between isotopic values of prey sources were inspected.Species (or groups of species) with large negative correlation values (> 0.5) indicate that multiple solutions exists with either one of the species, which is reflected by larger credible intervals in the resulting posterior distribution.Seasonal variations of trophic interactions and isotope values were considered by running separate models for each season.To better visualize the importance of mesozooplankton prey according to their size-class, contributions where then aggregated a posteriori as a function of the maximum size of the species pools and corresponding to Seston (or POM), " < 1 mm" for Nauplii Cirripedia, Euterpina acutifrons and Cypris; "1 to 1. www.nature.com/scientificreports/MixSIAR models were run using predator-prey seasonal TEF values (mean ± sd) previously calculated by IsoWeb (Supplementary S4).TEF variation across links was estimated, assuming that TEFs follow a normal distribution with a mean of 0.8 for carbon and 2.3 for nitrogen, as proposed for zooplankton food-webs 63,64 .Models were run seasonally with the following parameters: 106 chain length, 50 k burn-ins, and thin number 500 for three parallel Markov Chain Monte Carlo (MCMC) chains.Convergence was assessed using the Gelman-Rubin test (Gelman et al., 2014).MixSIAR models were run under the "very long" (10 6 chain length, burn = 500.000,thin = 500, chains = 3) or "extreme" setting for complex model with a large number of prey (3*10 6 chain length, burn = 1.500.000thin = 500, chains = 3).Convergence was assessed using the default MixSIAR diagnostic Gelman-Rubin and Geweke tests (see Supplementary S5 for the complete posterior distribution outputs).

Sources of variation of baseline SI values
Stable isotope values at the base of the planktonic food-web (i.e.seston) were highly variable and ranged from − 25.8‰ to − 16.2‰ for δ 13 C and from 3.4‰ to 9.5‰ for δ 15 N (Supplementary Fig. S1).Similarly, C:N ratios varied among seasons with lower values observed in spring (mean 9.13 ± 6.37) and higher values observed during autumn (12.4 ± 10.36) and winter (16.60 ± 12.53).GAM models for both nitrogen and carbon indicated that there was a non-linear relationship between the isotopic values and the spatial coordinates (p values < 0.001).In terms of the seasonal variability, both seasons, spring and winter were found to have significant negative effects on δ 15 N, indicating that δ 15 N tends to be lower in these seasons compared to autumn (p value < 0.01) (Table 2).Similarly, carbon values in winter were lower than in the other seasons.Further examination of the spatial patterns of estimated baseline values revealed that baseline estimates of δ 13 C are also significantly influenced by depth, while estimates of δ 15 N are explained by the interaction between depth and SPM, serving as a proxy for river influence (Supplementary S1).
The comparison of the isotopic niche space occupied by the plankton community (Fig. 3) highlighted different planktonic food-web architectures.In winter, the isotopic niche space is positioned to the right (higher δ 13 C adjusted ), reflecting an influence of benthic-derived carbon sources.Species with the highest isotopic divergence were cumaceans and E. acutifrons for zooplankton (lower end of δ 15 N adjusted ), and fish larvae of P. platessa and Ammodytidae for higher δ 15 N adjusted values.On the contrary, in spring, the isotopic niche space was centered around zero, thus reflecting a dominance of pelagic carbon-derived energy.Spring was also characterized by the higher range of δ 15 N adjusted values, with minimum values for the copepod E. acutifrons, nauplii Cirripedia, and cumaceans, and maximum values for fish larvae M. merlangius, T. luscus, and Gobiidae, indicating a higher vertically structured community (i.e., a greater number of trophic levels).The smallest isotopic niche space was found in autumn, characterized by a relatively large range of δ 13 C adjusted values, indicating the influence of both pelagic and benthic-derived carbon.Species displayed a more generalist strategy, as reflected by the lower range of δ 15 N adjusted values.

Table 2. Summary of Generalized Additive Models (GAMs) on baseline isotopic values. The table presents
smooth terms and parametric coefficients, along with estimated degrees of freedom (edf), reference degrees of freedom (Ref.df),F-values, and associated p values.For Nitrogen, the model explained 63.1% of the deviance with an R 2 of 0.58, while for Carbon, the model explained 52.6% of the deviance with an R 2 of 0.43.G¹³C adjusted

G¹⁵N adjusted
. Seasonal variation of the isotopic niche space occupied by the plankton community in the EEC.The black polygon illustrates the overall theoretical isotopic niche space, which can be compared to the seasonal realized niches (in blue, green, and orange for winter, spring, and autumn, respectively).Seasonal data are represented as black dots, and species at the edges, reflecting those with a higher trophic divergence, are identified.www.nature.com/scientificreports/

Food-web topology and energy fluxes
Size-structure seasonal variation of trophic pathways (MixSIAR models) Ichthyoplankton assemblages in winter were characterized by young-larval stages of herring (12 ± 3.4 SD mm SL) and plaice (8.85 ± 1.12 mm SL) and older stages of sardine (21.60 ± 1.95 mm SL).Among other carnivorous plankton, mysids and chaetognaths were also frequently encountered.Analysis of diet by size-classes showed that the main contributors for all winter species (representing between 51 to 82% of the diet) were between 1 to 1.5 mm in size corresponding to the copepods A. clausi, Paracalanus/Pseudocalanus spp, as well as D. anglicus.POM represented less than 10% of the diet for all species highlighting their predatory nature.Mysids diet was dominated by prey between 1 to 2 mm (both size classes accounting for ~ 66% of the diet) (Fig. 5).In spring, ichtyoplankton assemblages were characterized by older herring larvae (16 ± 3.8 mm SL).For herring, a clear pattern was found with a preference (44% ± 0.14 of the diet) for larger prey (> 2 mm) and a negligible contribution of POM that only accounted for 4% of the diet.A similar pattern was found for dab larvae (10 ± 2.6 mm with a greater contribution of larger prey (66% ± 9% of the diet) and only negligible inputs from POM (~ 2%).For sprat larvae (20 ± 3.5% mm SL), there was an increased contribution from 24% for small zooplankton (< 1 mm) to 36% for medium-size zooplankton between 1.5 to 2 mm.Species > 2 mm contributed to 13% of the diet.Callionymus spp.Larvae diet was dominated by prey between 1 to 2 mm and mysids had a preference for prey > 1.5 mm.In autumn, Sardine larvae of similar sizes as those in winter were found (~ 21 mm SL).Contributions by size classes showed that main prey (49% ± 14% of the diet) were around 1 to 1.5 mm length.A similar patterns was found in Chaetognaths (~ 7-8 mm TL) that presented a pattern dominated by prey of intermediate sizes.

Discussion
In this study, we investigated the variability in the planktonic food-web structure for a coastal ecosystem (the EEC), including mesozooplankton and fish larvae.We explored the isotopic niche space used by the plankton and its relationship with species' size.Additionally, we employed stable isotopes mixing models to explore how different food-web architectures transferred up to higher trophic levels.Our results showed that food-web structure varied seasonally with size and highlighted different feeding patterns (trophic redundancy vs. trophic  divergence) among seasons and species (Table 4).In the following sections, we discuss our results in light of (1) how the different structures and functions relate to changes in productivity and environmental drivers and ( 2) implications for the energy transfer to predatory plankton and in particular to fish larvae.Finally, (3) we discuss some of the remaining knowledge gaps and comment on what is needed to move forward.

Seasonal variation on planktonic food-webs
Following previous studies 19,21,65 , season emerged as a significant driver of variation for the mesozooplankton community structure.Overall, some of the variation could be attributed to seasonal differences in species composition contributing to the different size classes.However, a considerable number of species of similar sizes (~ 44%) were present throughout the year and exhibited substantial variation in their stable isotope composition, which indicates that the variability is also likely due to seasonal changes in species' diets.Flexibility of feeding strategies is a well-known and common pattern for copepods that are usually considered opportunistic omnivores.
As an example, in spring, the planktonic harpacticoid copepod E. acutifrons was characterized by the lowest (and negative) values of δ 15 N adjusted , indicating that their isotopic composition was lower than the seston (proxy of POM) used as a baseline.This could be explained by the capacity of E. acutifrons to display selective feeding in ecosystems or instances with high suspended particulate matter levels 66 (e.g., from increased rainfall or river discharge due to seasonal weather patterns) like those that characterize the EEC.On the contrary, nitrogen isotope values of E. acutifrons in winter and autumn were centered on zero, thus reflecting a more herbivorous diet when suspended particulate matter was low.Similar patterns were found for other copepods, such as T. longicornis (a suspension feeder), or A. clausi (an ambush feeder).These species are considered mostly herbivorous but they can prey (and even preferentially select) heterotrophic protists (ciliates) when phytoplankton concentrations are low 31,67,68 .As expected, the highest δ 15 N adjusted values were recorded in carnivore plankton (i.e., chaetognaths, mysids and fish larvae) for all seasons 21 , thus reflecting a higher trophic position for these species.However, the range of δ 15 N adjusted for the autumn food-web was the lowest among all seasons (as reflected by the smaller isotopic niche space), and values for fish larvae were of the same order as some other mesozooplankton species, such as the shrimp C. crangon or the copepod A. clausi.A large range of δ 13 C adjusted values was found in winter and autumn, indicating diversity in the origin of carbon sources (benthic and pelagic).Maximum values of δ 13 C adjusted were found in cumaceans that had the highest trophic divergence at all seasons.Cumaceans feed mainly on microorganisms and organic material from the sediment, thus reflecting benthic organic matter consumption 69 .The presence of meroplanktonic larvae of benthic species further contributes to interactions between the plankton and the benthos 70 .Our results highlight the significance of benthic-pelagic coupling, a crucial pattern recognized in coastal ecosystems, especially in relatively shallow and well-mixed waters like the Eastern English Channel 59,71 .While previous studies focused on fish, our findings extend these vital connections to planktonic organisms, underlining the pervasive influence of benthic-pelagic coupling across multiple trophic levels.Seasonal differences were also reflected in the size-structure of planktonic food-webs as observed by the different slopes and seasonal estimates of our LMEM models.Higher δ 15 N values with increasing size are frequently explained as reflecting size-related feeding patterns in marine plankton food-webs 23 .The relationship between nitrogen and size was strongest in spring, suggesting that the energy derived from phytoplankton blooms resulted in a more size-structured food-web and more specific predatory diets (higher trophic divergence as seen in the isotopic niche space).Similar trends have been found in tropical and subtropical regions where species tended towards more carnivorous feeding strategies leading to a higher vertical trophic structure (i.e., larger range of nitrogen isotope values or trophic levels), in periods of high Chl a (cold, non-stratified water) than during less productive seasons 24,72 .Additionally, strong trophic vertical structures during productive seasons have been hypothesized as the result of the accumulation of biomass and stronger microbial food-webs that increased food-chain length 19 .Conversely, the gentle slope of δ 15 N values with size in autumn, where consumers and prey displayed very similar values, suggests that the planktonic food-web is more dependent on recycled production 10 .These results are supported by the higher C:N ratios and more depleted δ 13 C values in autumn and winter, indicating a more processed material, possibly resulting from decomposition and recycling within the marine system and potentially from a benthic origin.Recent studies have shown that at least some protists may exhibit variations in 15 N trophic enrichment that deviate from the well-established patterns observed in metazoan consumers.Therefore, it is possible that when microbial activity dominates the energy pathways in the plankton, it www.nature.com/scientificreports/leads to lower δ 15 N values for consumers 26 .Conversely, lower C:N ratio and less depleted δ 13 C values observed in spring suggest a period of increased primary productivity, likely due to phytoplankton blooms.Linked to the isotopic niche space occupied in each season, trophic diversity was highest in spring and lowest in autumn (generalist species-trophic and functional redundancy).Our results are in agreement with previous studies on the Mediterranean suggesting that food overlap (or trophic redundancy) among zooplankton species and size classes seems higher during less productive (summer-autumn) than during the high productivity seasons (late winter-spring) 65 .In line with recent observations for the north Atlantic, our data highlights that zooplankton food-webs are organized in complex trophic structures that are not easily summarized into 2-3 functional groups disregarding seasonal inter-and intraspecies variation in feeding patterns.The seasonal variations in the planktonic food-web architecture from spring to winter and autumn align with previous observations of a continuum of trophic structures where 'herbivorous-based food-webs' vs 'microbial-based food-webs' represent only extreme configurations of the transient nature of a single planktonic food web.These configurations vary seasonally, depending on nutrient availability, phytoplankton production, and bacterial activity 73 (Fig. 6).

Energy transfer to predatory plankton
Fish larvae are key players in food-web and population dynamics.Their survival is considered one of the main processes influencing stock recruitment variability 1,2 and is closely to their capacity to capture suitable prey 3,74,75 .Considered mainly omnivorous/carnivorous feeders, research over the past decades has shown that fish larvae can display different feeding strategies in order to reduce potential interspecific competition with other predators 76,77 .Similarly, ontogenetic changes linked with body size (and by extension mouth size) suggest that most species switch to bigger prey as they grow 78 .However, prey selection based on species-size, and tradeoffs between prey availability and capturability are also part of fish larvae foraging strategies 32 .In this study, we investigated possible carry-over effects of the different seasonal planktonic food-webs structures on fish larvae and other carnivorous plankton.
In winter, the ichthyoplankton assemblage of the EEC is dominated by herring larvae after spawning of the Downs herring component.Plaice and sardine larvae are also frequently encountered, although to a lesser extent.The former species have been reported to be omnivorous at their larval stage based on stomach content analysis 30,79,80 , and variation of diet among regions has led researchers to believe that they might feed on the most abundant prey 81,82 .The use of stable isotope mixing models is complementary to stomach contents analyses, as the former reflects the assimilated diet while the latter informs on ingested diet.Our results agree with previous observations indicating that winter predators (including mysids and chaetognaths) can feed on a wide variety www.nature.com/scientificreports/ of prey sources.Theoretically, under food-limited conditions, fish larvae cannot afford to select prey and should ingest a wider range of prey sizes 83 .However, experimental studies show that at colder temperatures (as expected in winter during low production) larvae depend more heavily on optimal prey sizes 32 .This is in agreement with stable isotopes mixing models indicating that small copepods (~ 1 to 1.5 mm) dominate the diet of winter predators.Unfortunately, we had no data on zooplankton abundance by size-classes to test if small copepods correspond to the most abundant size-class, thus supporting that fish larvae in winter behave as opportunistic predators.Similarly, the contribution of phytoplankton to fish larvae has been reported as an important food for first-feeding and young larvae that can use diatoms as a type of initial or exploratory food source, to establish their feeding behavior 82,83 .In our study, seston (POM) appeared as a negligible contributor to the diet (< 10% of the diet) suggesting that older larvae are more carnivorous than omnivorous and that the nutritional value of POM is rather limited.The surprisingly large contribution (> 75%) of small copepods for plaice larvae might however, be an overestimation, as data on one of the main potential prey of plaice, the appendicularian Oikopleura dioica 81 was missing (see below discussion on missing species).However, previous studies have shown that small copepods such as A. clausi and Para-Pseudocalanus have poor escape capabilities 84 that might lead to a positive selection by fish larvae 78 .
Ichthyoplankton in spring was more diverse when compared to the other seasons, and so we expected some differences in feeding patterns among fish larvae to lower possible interspecific competition.Following the strong size-structure of the zooplankton food-web, dietary patterns in herring and dab larvae indicate a dominance of bigger prey in the diet (> 1.5 mm) and underline the importance of large copepods such as C. helgolandicus for the transfer of energy to higher trophic levels.Herring larvae in spring are bigger when compared to individuals collected in winter and are expected to be more successful at capturing bigger prey.Still, larger prey also nated the diet of the dab larvae (~ 10 mm).Contrary to previous gut content analysis that found that dab feeds mainly on small items (nauplii and copepodites of T. longicornis) 85 , our results suggest that small prey represent less than 10% of the assimilated diet in the EEC.Sprat and dragonet larvae did not seem to have a preferred or dominant size-class.Results for sprat are in agreement with previous studies that show that trophic niche breadth increased with larval' size from newly hatched to pre-schooling larvae (~ 16 mm) but then remained unchanged with sprat larvae feeding (and selecting) prey of different size-classes such as Acartia spp and C. hamatus 83 .The contribution of larger prey such as C. helgolandicus only represented ~ 13 to 24% of the diet, which suggest that some predator species might have adapted their feeding strategy to a more generalist diet to lower possible inter-specific competition and avoid trophic niche overlap.In agreement with previous studies, mysids appeared as carnivorous, feeding on similar zooplankton prey as fish larvae.The only exception was the negligible contribution of A. clausi (~ 5%, see Supplementary S5) that seems to be rejected by mysids even when the prey is abundant in the water 86,87 .
In autumn, larvae of several fish species were recorded but unfortunately only sardine larvae were collected in sufficient numbers for stable isotopes mixing models.Transfer of energy to predatory plankton was only explored for sardine and chaetognaths that were also frequently encountered.Both predators fed on a variety of species from different size-classes.Contrary to herring, sardine larvae have a late spring-summer and autumn spawning season 88 so that larvae collected at both periods were of similar sizes.Similar to the winter pattern, small prey (around 1 to 1.5 mm) dominated the diet of sardine larvae.In autumn, contrary to other seasons, δ 15 N adjusted values for chaetognaths and mysids were similar to those of fish larvae, suggesting that both groups share a similar trophic level and might feed on similar resources.These results concurs with the smaller isotopic niche space and higher trophic redundancy for the plankton food-web observed during autumn.
Overall, seasonal patterns of planktonic food-webs seem to propagate to upper trophic levels including fish, in particular at low productive seasons.A recent study 89 , looking into the plasticity of adult fish assemblages in the EEC during autumn and winter, also showed a reduction of the isotopic niche space and number of trophic levels in autumn, and a higher vertical structure in winter.The authors suggested that changes in feeding strategies were probably the result of differences on primary production (leading to changes in prey abundances and possible competition or niche overlap).Although we have no data on zooplankton abundance to confirm this observation, foraging on similar prey in autumn suggests that there is no clear density dependence, which leads to trophic similarity 90 .In spring, different feeding patterns related to prey size suggest that other factors than consumer' body size (e.g., resource partitioning or competition) influence larval feeding strategies during productive periods.

Remaining knowledge gaps and future directions
Other sources of variation: There are multiple sources of variation and uncertainties when using stable isotopes to elucidate trophic patterns [91][92][93] .For instance, possible inter-annual variations on plankton stable isotopes were not explored in this study because of data limitations.However, studies have shown that baseline values can vary spatially (isoscapes) but that these spatial patterns are stable from year to year (summer values over 10 years in the North Sea) 94 .Additionally, recent analysis of zooplankton community and size structure between 1991 to 2013 in the EEC in winter, showed that patterns (in terms of community composition, abundance and size-structure) where relatively stable over time within the region of our study 42 .Spatial and environmental patterns can also influence isotopic values at the base of the food-web.The distinct spatial trends in seston δ 13 C and δ 15 N are likely a result of a complex interplay of factors, including resource availability for phytoplankton, phytoplankton community structure, and the mixing of organic matter from various sources.These factors exhibit strong spatial gradients from the coast to offshore and/or from the distance to the river plume [95][96][97] .This clearly paves the way for further investigation into the factors influencing baseline reference values.Notably, our results align with isoscapes estimated using alternative modeling approaches (i.e.integrated nested Laplace approximation, INLA), as applied in a study of predatory gelatinous zooplankton in the EEC during 2015-2016 97 .
In that context, temperature, phytoplankton taxonomy, terrestrial nutrient input (at a relatively local scale), and mixing degree emerged as pivotal factors shaping isotopic baseline structures 98 .Taken together, these studies suggest that major changes in both zooplankton communities and stable isotopes in the EEC are likely driven by changes in temperature (probably an indirect link) and productivity, supporting our statement that season is the main driver of variability for planktonic food-webs at relatively low spatial scales.
Missing species: Copepods represent the majority (~ 90%) of the mesozooplankton in the EEC 42 .However, gelatinous zooplankton (Cnidaria, Ctenophora, Tunicata) are also frequently encountered 99,100 and can occasionally occur in large numbers with biomass exceeding that of fish in oligotrophic waters 21 .Unfortunately, samples of gelatinous zooplankton were not preserved from our surveys and so knowledge on their trophic dynamics and seasonal variability for the EEC planktonic food-web remains limited.Even though some species are too fragile and difficult to identify, samples of large cnidarians and ctenophores can help elucidate the trophic structure of gelatinous species and possible predation or niche overlap with fish larvae 101 .Joint information on crustaceans, gelatinous zooplankton and ichthyoplankton can be used as indicators of energy flow and trophic pathways, which should inform on how planktonic communities respond to environmental changes.Such information is required to inform several management descriptors (e.g.OSPAR indicators, Marine Strategy Directive-D1-Biological Diversity D4-Marine Food-webs, https:// oap.ospar.org/).

Conclusion
In conclusion, our study highlights significant seasonal variability in the planktonic food-web in the Eastern English Channel and Southern Bight of the North Sea.These dynamics, intricately linked to fluctuations of suspended particulate organic matter and primary production, give rise to distinct variations in isotopic niche spaces and trophic structure (both vertically in terms of trophic levels, and horizontally in terms of carbon source variability), as well as size-structural patterns.Remarkably, the food-web architecture of lower trophic levels propagate upwards to carnivorous plankton including fish larvae, and even higher up to adult fish.This emphasizes the pivotal role of bottom-up control in shaping coastal systems like the EEC, and the need for regular, long-term monitoring of lower trophic levels across the full seasonal and spatial gradients of any given management area.A better understanding on how environmental parameters shape trophic transfers is essential to predict how planktonic food-webs will respond to global change scenarios.Furthermore, the intrinsic trophic link between mesoplankton and fish larvae (alongside other vital resources for sustaining human consumption) underscores the necessity for enhanced management and an ecosystem-based approach that includes planktonic species based on life-history traits and size spectra.

Figure 1 .
Figure 1.Study area including the Eastern English Channel and Southern Bight of the North Sea.Sampling stations in winter (blue), spring (green) and autumn (yellow) and main rivers along the French coast are indicated.

Figure 2 .
Figure 2. Variability (seasonal, inter, and intra-specific) of plankton isotopes (left, δ 13 C; right δ 15 N; top, mesozooplankton; bottom, fish larvae) in the EEC.Species are ordered based on their averaged δ 15 Nadjusted values.Seasonal mean values are illustrated by dots.Unique values (if n = 1) are illustrated by squares.

Figure 4 .
Figure 4. Seasonal and size effects on δ15 Nadjusted and δ13 Cadjusted values of the plankton community.Lines illustrate predictions from the Linear Mixed-Effects Models (LMEM) with 'Season' and 'Size' as fixed factors.The best carbon model includes 'Species' as a random effect, while the nitrogen model incorporates variation in size within species ('log_size | Species').Zooplankton values are illustrated as triangles, and fish larvae are represented as circles, colored according to the season.

Figure 5 .
Figure 5. Diet composition (% diet) of plankton size-classes to predatory plankton.Values represent mean values and standard deviation posterior distributions of the MixSIAR models.Smooth dashed lines are for illustration purposes only and highlight main patterns or dominant size-classes to the diet.

Figure 6 .
Figure 6.Schematic representation of the continuum of trophic structures for the EEC planktonic food-web.

Table 1 .
15mber of measurements by species and by season (pools of 1 to 100 individuals).Mean and standard deviation (sd) of size (total body length, mm) and stable isotopes values for δ15N and δ 13 C.A total of 552 measurements of planktonic species are recorded.

Table 3 .
Linear Mixed Effect Models (LMEM) examining the influence of season, size, and species on zooplankton values.Marginal R 2 and conditional R 2 are used to assess the proportion of variance explained by fixed and random effects.σ 2 represents unexplained variability in the response variable, while τ00 accounts for variability between species levels.p values are based on conditional F-tests using the Kenward-Roger approximation for degrees of freedom and the pbkrtest-R package.

Table 4 .
Summary of main findings, highlighting different food-web architecture indicators and their seasonal variation.